Non-Gaussian Fluctuations in Biased Resistor Networks: Size Effects 

versus Universal Behavior 

C. Pennetta 1 , E.Alfinito 1 , L. Reggiani 1 , S. Ruffo 2 

1 Dipartimento di Ingegneria dell'Innovazione, Universitd di Lecce and 
National Nanotechnology Laboratory-INFM, Via Arnesano, 73100 Lecce, Italy. 

2 CSDC, INFN and Dipartimento di Energetica "Sergio Stecco", 
Universitd di Firenze, Via S. Marta, 3, Firenze, 50139, Italy. 

keywords: Non-Gaussian Fluctuations, Resistor networks, Nonequilibrium stationary states, 
Disordered materials 

PACS: 05.40-a, 05.70.Ln, 64.60.Fr, 72.80 Ng 

Abstract 

We study the distribution of the resistance fluctuations of biased resistor networks in 
nonequilibrium steady states. The stationary conditions arise from the competition be- 
tween two stochastic and biased processes of breaking and recovery of the elementary resis- 
tors. The fluctuations of the network resistance are calculated by Monte Carlo simulations 
which are performed for different values of the applied current, for networks of different size 
and shape and by considering different levels of intrinsic disorder. The distribution of the 
resistance fluctuations generally exhibits relevant deviations from Gaussianity, in particular 
when the current approaches the threshold of electrical breakdown. For two-dimensional 
systems we have shown that this non-Gaussianity is in general related to finite size effects, 
thus it vanishes in the thermodynamic limit, with the remarkable exception of highly dis- 
ordered networks. For these systems, close to the critical point of the conductor- insulator 
transition, non-Gaussianity persists in the large size limit and it is well described by the 
universal Bramwell-Holdsworth-Pinton distribution. In particular, here we analyze the role 
of the shape of the network on the distribution of the resistance fluctuations. Precisely, we 
consider quasi-one-dimensional networks elongated along the direction of the applied current 
or trasversal to it. A significant anisotropy is found for the properties of the distribution. 
These results apply to conducting thin films or wires with granular structure stressed by 
high current densities. 

1 Introduction and Model 

Strongly correlated systems usually exhibit non-Gaussian distributions of the fluctuations of global 
quantities, as a consequence of the violation of the validity conditions of the central-limit theorem. 
Since correlations become important near the critical points of phase transitions, non-Gaussian 
fluctuations are usually observed near criticality [1, 2, 3, 4, 5, 6, 7, 8, 9]. In these conditions, the 
self-similarity of the system over all the scales, from a characteristic microscopic length up to the 
size of the system, has important implications on the fluctuation distribution [2, 3, 4, 5, 6, 7, 8, 9]. 
On the other hand, far from criticality, the correlations among different elements of the systems 
can also be important. This is particularly true for systems in non-equilibrium stationary states, 
where non-Gaussian fluctuations are frequently present [1, 10, 11, 12, 13, 14, 15]. Therefore 
the study of non-Gaussian fluctuations and of their link with other features of the system can 
provide new insights into basic properties of complex systems [2, 3, 4, 5, 6, 7, 8, 9, 11]. On this 
respect, the observation made few years ago by Bramwell, Holdsworth and Pinton (BHP) [2] of 
a common behavior of the distribution of the fluctuations in two very different systems, (the 



power-consumption fluctuations in confined turbulent-flow experiments and the magnetization 
fluctuations in the two-dimensional XY model in the spin- wave regime at low temperature [2, 
3, 4]), has given rise to several intriguing questions about the origin of this common behavior, 
stimulating many other experimental, analytical and numerical studies. Successive findings have 
highlighted that many scale invariant systems display the same functional form for the distribution 
of the fluctuations [3, 7, 5, 13, 14, 15]. Very recently, new light on these puzzling observations 
has been given by Clusel et al. [16]. These authors, on the basis of a study of the fluctuation 
properties of the 2D XY model, have proposed a criterion for universality-class-independent 
critical fluctuations [16]. Actually, in the relatively simple case of the 2D XY model it is possible 
a complete understanding of fluctuation phenomena. This is not possible for nonequilibrium 
systems due to the lack of microscopic theories. Thus, for these systems, we can rely only on 
phenomenological observations and on analogies with better understood systems. 

Here, we study the distribution of the resistance fluctuations of biased resistor networks in 
nonequilibrium stationary states [17]. Networks of different size and shape and with different levels 
of internal disorder are considered. The resistance fluctuations are calculated by Monte Carlo 
simulations for currents close to the threshold for electrical breakdown. This last phenomenon 
consists of an irreversible increase of the resistance, occurring in conducting materials stressed by 
high current densities and it is associated with a conductor-insulator transition [18, 19, 20, 21, 
22, 23]. In our study we make use of the Stationary and Biased Resistor Network (SBRN) model 
[24, 25]. This model provides a good description of many features associated with the electrical 
instability of composite materials [20, 22, 24] and with the electromigration damage of metallic 
lines [17, 21], two important classes of breakdown phenomena. 

We describe a thin conducting film with granular structure of length L, width W and thickness 
t h <C W, L as a 2D resistor network of rectangular shape and square-lattice structure [17]. The 
network of resistance R is made by Nl and Nyy resistors in the length and width directions 
respectively. Thus, the total number of resistors in the network (excluding the contacts) is: 
N tot = 2N L N W + N L — N w . The external bias (here a constant current /), is applied to the 
network through electrical contacts realized by perfectly conducting bars at the left and right 
hand sides of the network. The network lies on an insulating substrate at temperature T , acting 
as a thermal bath. Each resistor has two allowed states [21, 26]: (i) regular, corresponding to 
a resistance r reg ^ n (T n ) = r re f[l + a(T n — T re f)] and (ii) broken, corresponding to a resistance 
r Q p = 10 9 r re9 n (T ) = 10 9 r (resistors in this state will be called defects). In the above expression 
a is the temperature coefficient of the resistance (TCR), r re f and T re f are the reference values for 
the TCR and T n is the local temperature. The existence of temperature gradients due to current 
crowding and Joule heating effects is accounted for by taking the local temperature of the n-th 
resistor given by the following expression [21]: 

T n = T + A[r n i 2 n + (3/AN netg ) £ (r,i? - r n i 2 n )} (1) 

i=i 

where, i n is the current flowing in the nth resistor and N neig the number of its nearest neighbors 
over which the summation is performed. The parameter A represents the thermal resistance of 
each resistor and sets the importance of Joule heating effects. By taking the above expression 
for T n we are assuming an instantaneous thermalization of each resistor at the value T n [21, 26]. 
In the initial state of the network (no external bias) we take all the resistors identical (perfect 
network). We assume that two competing biased processes act to determine the evolution of 
the network [24, 25]. These two processes consist of stochastic transitions between the two 
possible states of each resistor and they occur with thermally activated probabilities [26] : Won = 
exp[—E D /kBT n ] and W Rn = exp[—E R /kBT n ], characterized by the two energies, E D and E R 




Figure 1: Pattern of a network 12 x 50 stressed by a current density j = 0.32 mA. The grey 
boxes show the backbone of the network, the black ones the "dangling bonds" (branches with 
zero-current, while the missing boxes correspond to the broken resistors. This pattern has been 
calculated at t — 4 x 10 4 (time expressed in iteration steps). 



(where ks is the Boltzmann constant). The time evolution of the network is obtained by Monte 
Carlo simulations which update the network resistance after breaking and recovery processes, 
according to an iterative procedure described in detail in Ref. [24]. The sequence of successive 
configurations provides a resistance signal, R(t), after an appropriate calibration of the time scale. 
Depending on the stress conditions (/ and To) and on the network parameters (size, activation 
energies and other parameters dependent on the material, like r re f, a and A), the network either 
reaches a stationary state or undergoes an irreversible electrical failure [17, 24, 25]. This latter 
possibility is associated with the achievement of the percolation threshold, p c , for the fraction 
of broken resistors. Therefore, for a given network at a given temperature, a threshold current 
value, Ib, exists above which electrical breakdown occurs [24]. For values of the current below 
this threshold, the steady state of the network is characterized by fluctuations of the fraction of 
broken resistors, Sp, and of the resistance, 5R, around their respective average values < p > and 
< R>. In particular, we underline that in the vanishing current limit (random percolation) [27], 
the ratio A = (E D — E R )/k B T determines the average fraction of defects and thus the level of 
intrinsic disorder inside the network [27]. In the following we analyze the results of simulations 
performed by considering networks of different size and shape stressed at room temperature, 
T = 300 K, by a current density j = I/N w = 0.32 mA. We have taken: a = 3.6 x 10~ 3 K -1 , 
T re/ = 273 K, r ref = 0.048 Q, A = 2.7 x 10 8 K/W, E D = 0.41 eV and E R = 0.35 eV. This choice 
of the parameters is appropriate to describe the behavior under electromigration of metallic lines 
of Al-0.5%Cu studied in Ref. [17] and it corresponds to studying a network with an intermediate 
level of intrinsic disorder. 



2 RESULTS AND CONCLUSIONS 

Figure 1 displays the pattern of a network 12 x 50 calculated at a given time, t — 4 x 10 4 , 
(expressed in iteration steps) in the stationary regime of the network, i.e. for t > T re i, where 
r re i ps 8 x 10 3 is the relaxation time for the achievement of the nonequilibrium stationary state. 
The network in this figure is stressed by a current density (j = 0.32 mA) close to the breakdown 
value, jb = Ib/Nw- The resistance evolution for the same network is reported in Fig. 2. In this 
figure the grey line shows the average value of the resistance, < R >. We note that both the 
average resistance and the relative variance of the resistance fluctuations, < {SR) 2 > / < R > 2 , 
depend on j. A detailed analysis of the behavior of these two quantities as a function of the 
current can be found in Refs. [24, 25]. In previous works [13, 14, 15] we have analyzed the 
effects on the distribution of the resistance fluctuations of the biasing current [13], of the intrinsic 
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Figure 2: Resistance evolution of the network in Fig. 1. The time is expressed in arbitrary units 
(iteration steps), the resistance in Ohm. The grey line shows the average value of the resistance. 



disorder and of the size of the network [14, 15], by limiting ourself to discuss square networks. 
Here, we focus our discussion on shape effects: precisely we analyze the effect on the distribution 
of the fluctuations of scaling the size of the network separately in the two directions, i.e. of 
scaling separately the width, N w and the length, N L , of the network. Figure 3(a) shows the 
distributions of the resistance fluctuations obtained for two networks of size 12 x 50 (big circles) 
and 50 x 12 (triangles) stressed by the same current density. In this figure (and in the followings) 
we denote with <3> the probability density function (PDF) of the distribution and with a the root 
mean square deviation from the average value. This normalized representation, by making the 
distribution independent of its first and second moments, is particularly convenient to explore 
the functional form of a distribution [3]. A lin-log scale is adopted for convenience to plot the 
product a& as a function of (< R > —R)/a. The PDFs in Fig. 3 and all the others in this 
paper have been calculated by considering time series containing about 10 6 resistance values. For 
comparison, in Fig. 3 we also report the Gaussian distribution (dashed curve) and the BHP 
distribution (continuous curve) [2, 3]. The PDF obtained for the network 12 x 50 (corresponding 
to the signal in Fig. 2) exhibits a strong non-Gaussianity, well described by the BHP curve. By 
contrast, the PDF obtained for the network 50 x 12 is nearly Gaussian. At a first insight, this 
result can seem surprising: in fact the two networks are composed by nearly the same number of 
resistors, moreover the dissipated electric power per unit volume, RI 2 /(LW) cx j 2 , is the same in 
both cases. As a consequence, the average fraction of defects p ~ 0.19 is also the same. However, 
the percolation threshold pc is different for the two networks [17]. Therefore, the nearly Gaussian 
distribution of the 50 x 12 network is due to the higher value of p c (and thus to the higher stability) 
of this network [17]. For comparison, we report in Fig. 3(b) the PDFs calculated for two square 
networks 12 x 12 and 50 x 50 biased by the same current density. Again, the dissipated power 
density and the average fraction of defects are the same for both networks. However, for square 
N x N networks the percolation threshold is roughly independent of the size, even for biased 
percolation [17]. Thus, for these networks the higher instability and the stronger non-Gaussianity 
for decreasing N is mainly related with the increase in magnitude of the fluctuations associated 
with the smaller size [13, 14]. 

The normalized PDFs of the resistance fluctuations calculated for several networks of different 
size are reported in Fig. 4. Figure 4(a) displays PDFs obtained for networks elongated along 
the direction of the applied current (precisely networks of a given width, N w = 12, and with 




Figure 3: Normalized PDF of the resistance fluctuations of networks of different size and stressed 
by the same current density j = 0.32 mA. Precisely, in (a) the size is: 12 x 50 (big circles) and 
50 x 12 (triangles); in (b): 12 x 12 (stars), 50 x 50 (plus). The solid and dashed curves refer to 
the BHP and Gaussian distributions, respectively. 



increasing length, Nl = 50 4- 400), while Fig. 4(b) shows PDFs obtained for networks elongated 
in a direction trasversal to the current applied (precisely networks of a given length, Nl = 12, 
and with increasing width, Nw = 50 -j- 400. We can see that for trasversal networks the PDF 
is rather insensitive to the width and the small non-Gaussianity for small widths completely 
vanishes already for networks with N\y = 200. By contrast, for longitudinal networks, the PDF 
is sensitive to the length. However, it should be noted that the PDF obtained for Nl = 200 
practically overlaps with that obtained for Nl = 400 and both exhibit non-Gaussian tails. Since 
the correlation length, £, for these networks is estimated to be £ < 5, networks with N L = 400 
can be considered as infinitely long. Thus, Fig. 4(a) suggests a persistent non-Gaussianity for 
longitudinal networks in the limit Nl — > oo, associated with the finite size of the network in 
the transversal direction. Furthermore, the magnitude of this non-Gaussianity is expeceted to be 
controlled by the level of intrinsic disorder. 

In conclusions, we have studied the distribution of the resistance fluctuations of biased resistor 
networks in nonequilibrium stationary states. We have considered networks biased by currents 
close to the threshold of electrical breakdown. As a general trend, the distribution of the fluctua- 
tions is found to exhibit relevant deviations from Gaussianity, which are in general related to finite 
size effects [13, 14]. However, for systems close to the critical point of the conductor-insulator 
transition, the non-Gaussianity persists in the large size limit [13, 14] and it is well described 
by the universal Bramwell-Holdsworth-Pinton distribution. Furthermore, we have analyzed the 
role of the shape of the network on the distribution of the resistance fluctuations, by considering 
quasi-one-dimensional networks elongated along the direction of the applied current or trasversal 
to it. A significant anisotropy is found for the properties of the distribution. These results apply 
to conducting thin films or wires with granular structure stressed by high currents. 
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Figure 4: Normalized PDF of the resistance fluctuations of networks of different size and stressed 
by the same current density j = 0.32 mA. Precisely, in (a) the size is: 12 x 50 (big circles), 
12 x 100 (crosses), 12 x 200 (full squares) and 12 x 400 (down triangles); in (b) the size is: 50 x 12 
(triangles), 100 x 12 (small circles), 200 x 12 (full diamonds) and 400 x 12 (squares). The solid 
and dashed curves have the same meaning of Fig. 3. 
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